Heralded Polynomial-Time Quantum State Tomography 
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We describe an algorithm for quantum state tomography that converges in polynomial time to 
an estimate, together with a rigorous error bound on the fidelity between the estimate and the true 
state. The result suggests that state tomography on large quantum systems may be much more 
feasible than the exponential size of state space suggests. In many situations, the correctness of the 
state estimate can be certified from the data alone, with no a priori assumptions on the form of the 
measured state. 



Quantum state tomography (QST) — the process of 
estimating a quantum state using local measurements 
given a large number of copies — is an exceptionally 
daunting task. As it is typically formulated [1], simply to 
output an estimate for a generic state of n qubits would 
take exponential time in n, given that there are an ex- 
ponential number of coefficients in a generic state's de- 
scription. This is but one of several inefficiencies. Most 
quantum states have exponentially small amplitudes in 
almost every basis, so to distinguish any one of these 
amplitudes from zero takes exponential time. Assuming 
one were able to collect all of that data from an infor- 
mationally complete measurement, one is left with the 
intractable computational task of inverting the measured 
frequencies to find an estimate of the state. 

These barriers to performing efficient quantum tomog- 
raphy for generic states are not necessarily applicable to 
real-world situations, however, as the traditional repre- 
sentation of quantum states in terms of Hilbert space 
vectors and density matrices is in a sense too general. In- 
deed, states which occur in many practical situations are 
specified by a small number of parameters. An efficient 
description could be a polynomial-sized quantum circuit 
which outputs the state; or, in the case of thermody- 
namical equilibrium states, a bounded-interaction-rangc 
Hamiltonian and a temperature. This insight is not new: 
researchers in many-body physics and quantum informa- 
tion theory have found many classes of states which are 
described by a number of parameters scaling polynomi- 
ally in n and which closely approximate the kind of states 
found in particular physical systems [2-8]. 

However, the question of whether these restricted 
classes can be put to use in a tomographic setting has re- 
mained largely open. Any result closing this gap should 
address three points: 1. Choose a few-parameter set of 
states adapted to the physical nature of the system under 
consideration. 2. Specify a protocol which efficiently [22] 
identifies these parameters from a small and simple-to- 
implement set of measurements. 3. If possible, give an 
independent way of verifying that the protocol produced 
a faithful estimate of the physical state (without having 
to assume that the set chosen in the first step correctly 



describes the system). 

We briefly explain the third point. An exprimentalist 
should be wary of any tomographic procedure for which 
the success is predicated on some abstract technical con- 
dition (e.g., the system is well-described by a "matrix 
product state of bond dimension D"). Therefore, it is 
desirable that the protocol infers from the data alone 
whether or not it can guarantee a faithful reconstruction, 
and aborts if it cannot. We say that such a procedure 
heralds its own success. 

In this work, we address all three problems above by 
presenting an algorithm for heralded quantum state to- 
mography using a polynomial amount (poly(n)) of mea- 
surement and classical computation. 

The physical system we have in mind is one where the 
constituents are arranged in a one-dimensional configu- 
ration (e.g. ions in a linear trap [9]). It is highly plau- 
sible that in such a setting, correlations between neigh- 
boring qubits are much more pronounced (due to direct 
interaction) than correlations between distant systems 
(mediated e.g. by global fluctuations of control fields). 
A polynomially-sized class of states anticipating exactly 
this behavior has long been studied under the names of 
finitely correlated states (FCS) or matrix product states 
(MPS) [2, 10]. Below, we present the heralded MPS to- 
mography (H-MPS) algorithm. Given experimental data, 
it outputs either a polynomial-sized pure matrix product 
state together with a rigorous bound on the fidelity be- 
tween the true density matrix and the estimate, or else 
it returns a failure flag. This happens in particular when 
such an efficient description does not exist. At least in 
an idealized scenario, one can guarantee that an efficient 
MPS description will be found in polynomial time, given 
that it exists. 

In later sections, we sketch generalizations to higher 
dimensions, to highly mixed systems, and to quantum 
channels describing correlated noise. We emphasize that 
the H-MPS algorithm can be employed whenever one ex- 
pects that an MPS description for a system may be ap- 
propriate, whether or not it possesses a one-dimensional 
geometry. Examples of MPS include the GHZ, W, clus- 
ter, and AKLT states. Our results also address an open 
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question of Aaronson [11] by demonstrating that MPS 
can be learned efficiently. 

To find an unknown quantum state p of n qudits (d- 
dimensional quantum systems), we propose the following 
procedure. One selects a number k <C n, such that ob- 
taining an estimate for the state of 2k contiguous qudits 
using standard tomographic techniques [1] is feasible. We 
denote the reduction of p to sites + + (2k — 1) 

by pj , and we use <jj to refer to a tomographic estimate 
of this state. Assume that by standard methods one ob- 
tains a confidence interval such that, for all j, the trace 
norm satisfies 

\\Pj - CTjIltr < €j (1) 

with high confidence. 

The data {o~j,6j} are then passed to the H-MPS al- 
gorithm. It either declares a failure and aborts, or else 
outputs an efficiently dcscribable state |W) and a number 
t such that the fidelity provably obeys 

F(\*),p)= ^MpW) > VT^ (2) 

whenever (1) holds. 

It may come as a surprise that such a certificate is 
possible. After all, there are exponentially many degrees 
of freedom which remain unmeasured in the procedure. 
How can we guarantee that there are no two orthogonal 
states which happen to coincide on the tiny number of 
coefficients we have obtained information about? Geo- 
metrically, the reason is that a generic MPS is "locally 
exposed" [2] , i.e. it lies at an extreme point of the space 
of quantum states, whose tangent plane corresponds to a 
local functional. In physical terms: with every MPS, one 
can associate a gapped, local Hamiltonian, which acts as 
a witness for that state [2, 10]. The main technical is- 
sue solved below is finding a simple way of bounding this 
gap for finite, non-translationally invariant and possibly 
degenerate systems. 

Before describing the algorithm, we pause to recall the 
definition of MPS [2, 10]. A state vector |*) on n qudits 
is an MPS with bond- dimension D if there are nd complex 
Dx D matrices {Aj }, s = 1, . . . , d; j = 1, . . . , n such that 

d 

|*) = Yl Tr[A^A?---A^]\ Sl , S2 ,..., Sn ). (3) 

si ,...,S Tl — 1 

Here, D is a free parameter which describes the complex- 
ity of the model. It is well-known that any state has an 
MPS description if one allows D to scale exponentially 
with n. However, the motivation for the seemingly ad- 
hoc definition (3) stems from the fact that many natural 
states are well-approximate by an MPS with small D. 
Note there are only nD 2 parameters specifying 

H-MPS. The algorithm consists of two steps. In the 
first stage, we compute an MPS estimate from the exper- 
imental data. The task of the second step is to certify the 



fidelity of the estimate. Because we will retrospectively 
verify the output of step one, it is not necessary to in- 
sist on theoretical guarantees on the performance of this 
part. Any ansatz, even an "educated guess", can lead 
to a rigorous result, as long as it will be vindicated by 
passing the second stage of H-MPS. Therefore, we give a 
list of algorithms for finding MPS estimates below, with- 
out definitely endorsing one over the others. A more 
detailed study of practical performances will be provided 
elsewhere. For clarity of presentation, we will restrict 
attention to "generic" states at first, treating singular 
situations in the end. 

The Variational Method. A natural way of obtaining 
FCS estimates from data is by way of a variational algo- 
rithm: note that the local reductions of an MPS can be 
computed efficiently. Therefore, it is feasible to change 
the matrices one at a time, so as to minimize the sum 
of the trace-norm differences between the measured es- 
timates and the computed local reductions. Since the 
target function is bounded from below and reduced in 
every step, the algorithm will certainly converge. The 
advantage of this method is that all the available infor- 
mation is utilized, and that it is inherently stable against 
small noise. On the negative side, it may be prone to 
run into local minima. (We emphasize, however, that 
the algorithm should not be judged by the high stan- 
dards of the variational methods employed in many-body 
physics, which reliably find the global minimum of en- 
ergy functions on tens of thousands of spins. A density 
operator even on ten three-level "qutrits" depends on 3.4 
billion coefficients. If the trace-norm distance can be re- 
liably minimized for such comparatively tiny systems, it 
would make tomography possible in otherwise intractable 
regimes - even if the algorithms fails to scale to the or- 
ders of magnitude which are standard today in numerical 
condensed matter physics). 

The Quantum 2- SAT Method. Recent developments in 
the context of the "quantum satisfiability problem" show 
that it is possible to directly and efficiently find the MPS 
description of the ground state space of a one-dimensional 
frustration-free Hamiltonian, provided that the following 
condition holds: for every 1 < j < n, the degeneracy of 
the restriction of the Hamiltonian to sites 1 to j does not 
exceed a uniform bound. A qubit version described in [12] 
can be generalized to arbitrary situations [13]. It easily 
follows that - at least in the idealized scenario of vanish- 
ing experimental errors - reconstructing an MPS descrip- 
tion given only information about the support of the re- 
duced density operators (on a sufficiently large neighbor- 
hood) is provably efficiently possible. While it is unclear 
whether this approach is well-suited for practical pur- 
poses, it shows that there are no computationally hard 
instances which would make any attempt of finding a 
general solution futile from the outset. This may be a 
surprising fact, given that hard instances are known to 
occur for more general frustration-free Hamiltonians with 
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unique MPS ground-states [14]. 

The DMRG Method. An MPS approximation with 
bond dimension D is only possible when the reductions pj 
are mainly supported on a space of dimension not larger 
than D 2 . Denote the projection operator on the space of 
the lowest {d 2k — D 2 ) eigenvalues of Oj by hj. One can 
use established methods like DMRG to find the ground 
state of this "empirical parent Hamiltonian" (see below) . 
Under suitable regularity conditions and assuming low 
noise levels, this procedure is guaranteed to yield good 
estimates, as will become clear later. Advantage: algo- 
rithms like DMRG perform extremely well under prac- 
tical conditions. Experience shows that they are fairly 
immune against local minima. The disadvantage is that 
only some parts of the measured information enters the 
procedure. 

The Educated Guess. The scenario most often encoun- 
tered in current experiments [9] is verifying the quality 
of a state preparation. Since we present a way of verify- 
ing correctness of a proposed solution, it makes sense to 
test a candidate MPS which did not come from measured 
data, but rather from the knowledge of which state one 
intended to create. This version of the problem would be 
more aptly described as state certification: simply use 
the MPS description of the target state to test the cer- 
tificate and get a lower bound on the fidelity of the state 
preparation. 

We proceed to the second stage. (The presentation will 
necessarily become somewhat more technical in the next 
paragraphs.) The estimate obtained in the first step is 
assumed to be given in MPS form (3). For simplicity of 
language, we group blocks of k contiguous qudits into a 
single site. In other words, we assume without loss of 
generality that k = 1 and that the individual systems 
are d fe -dimensional. (We will always assume that d k > 
D, since otherwise there are more unknown parameters 
in (3) than were measured.) To avoid numerically ill- 
conditioned inputs, we bring the MPS estimate into its 
standard form fulfilling J2 S ^ji^jV = 1> which is always 
possible by a sequence of singular value decompositions 
(SVD) [2, 10]. Next, compute the SVD of the maps 

d 

T 3 :X^ Tr[XA s /A s /£)\ Sj ,s j+1 ). (4) 

s j ■ s j+i= 1 

Generically, the singular values of the Tj's will be 
bounded away from zero [2, 10]. Should that not be 
the case, we declare an error and abort (this situation 
is discussed in more detail below) . If the estimate passed 
this first test, we are ready to construct and analyze its 
"parent Hamiltonian" [2, 10], which will serve as a wit- 
ness operator, certifying closeness of the estimate to the 
physical state. To that end, for every j <E {1, . . . , n — 1}, 
compute the reduced density matrix of the estimate \*S?) 
on sites j, Let hj be the projection operator onto the 
kernel of this density matrix. For every j £ {1, . . . , n — 1} 



compute the eigenvalues of hjhj + \hj. Let jj be the 
square-root of the largest eigenvalue not equal to one. 
Set 7 = maxj'Yj. The H-MPS algorithm concludes by 
outputting \^>) and t = J2j(T r (hj°~j) + e j)- 

In the following paragraph we prove the validity of 
Eq. (2). The technical hurdle to overcome is to lower- 
bound the gap in the spectrum of H = ^™ =1 hj. For 
non-degenerate translationally invariant (TI) systems, a 
(fairly involved) proof was given in [2]. Partly building 
on their ideas, we present a short argument below, which 
applies to the relevant case of finite, non-TI quantum 
states. 

By construction, H > 0, Hip = 0, so that the low- 
est eigenvalue of H is zero. We claim that the second- 
smallest eigenvalue gap(ff) is lower-bounded by (1 — 27). 
To prove the assertion, we start by following [2] and note 
that gap(£f) = max{A | H 2 > XH}. Compute: 

n-l n-2 

r2 = 12 h2 3 + 12^+ i+h ^ ih ^ + 12 
j=i j=i ij-j'i>i 

n-2 

> H + Y,[hj,h j+1 ] + , (5) 

where we used the fact that the hj are projections and 
omitted non-negative summands. Next, we will establish 

[hj,h j+1 ] + > ~lj{hj + h j+1 ). (6) 

We borrow some facts from basic Hilbert space theory 
[15, 16]: any two projection operators hj,hj+i can be 
brought simultaneously into block diagonal form, where 
on the ith block, hj and hj+i act like one-dimensional 
projections whose ranges are vectors enclosing the canon- 
ical angle One verifies that cos 2 9i is the ith eigenvalue 
of hjhj + \hj. Certainly (6) holds on the entire space if 
and only if that relation is true on every block. Thus, we 
have reduced the problem to the case of one-dimensional 
projectors in a two-dimensional space. Here, we can solve 
it by elementary means and find that (6) holds for non- 
trivial angles as long as 7 > cos 6i , which is true by con- 
struction. Plugging (6) into (5) gives 

H 2 > H -lYJ^i + h i+^ - _ 2 ^) ff ' 
3 

completing the proof of the gap estimate. To get the 
fidelity bound (2), label the eigenvectors of if by \Eq) = 
\^f),\E\),...\Edn—i). We know that the eigenvalue Eq 
equals zero. Since invcrtibility of the Tj's guarantees 
uniqueness of the ground-state [2, 10], the bound on the 
gap ensures that all other eigenvalues fulfill JSj > (1 — 27). 
Thus 

(Tr[V,] + ej) > Tr[Hp] 

3 

= Y,Ei(&\p\Ei) > (1 " 2 7)(1 - (*|p|*». 
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AKLT tomography, n=1 5 spins 
Fidelity density bound vs. total statistical error 




FIG. 1: Tomography of an AKLT state, a particular state of n 
spin-1 particles with an efficient matrix product description. 
The reconstruction method uses the DMRG method discussed 
in the main text for nearest-neighbor tomography only, and 
computes the error parameter r by using the DMRG estimate 
of the gap. The x-axis is the total error, E = . ej, and the 

y-axis is the lower bound on the fidelity density F 1 ^ 71 . The 
true fidelity of the recovered matrix product description was 
near 1 in all cases, suggesting that the method for finding the 
state performs much better than what the bound allows us to 
certify. The bound can be significantly improved if next-to- 
nearest-neighbor correlations are included. 

Finally, we remark that based on the ideas in [2], it can 
be shown that 7 — > exponentially fast as k is increased. 

Degeneracy and singular cases. While the certifica- 
tion algorithm above works for "most" MPS states, it 
may fail in certain singular cases. The simplest example 
is given by the family of GHZ-type states ^=(|0, . . . , 0) + 

e^|l, . . . , 1)). Since any local reduced density matrix is 
independent of <f>, it is impossible to distinguish the mem- 
bers of that family based on local information alone. (In- 
deed, for these states the -operators defined above will 
not be invertible - so that an error would have been de- 
clared). However, one may verify that the simple "string 
operator" a x ® . . . ® a x has expectation value cos <fi for 
the states above. Thus, it is reasonable to expect that PI- 
MPS can be extended, enabling it to learn (and certify) 
any MPS with bounded bond dimension. 

It turns out that this possible, as we will briefly sketch 
now. Indeed, if the algorithm described above reports an 
error (i.e. if the T^s fail to be invertible), one would pro- 
ceed as follows. Employing the Quantum 2-SAT Method, 
compute an MPS representation of the projector onto the 
ground space of the parent Hamiltonian of the estimate. 
It may happen that the dimension of this space grows ex- 
ponentially as more sites are taken into account. In this 
case, H-MPS will have definitely failed (in the limit of 
vanishing noise, this only happens if no MPS description 
with given bond dimension exists). However, in many 
relevant situations, the dimension will saturate at a fi- 
nite value. Re-examining the gap estimate above shows 
that (unlike the original proof in [2]), it remains valid 



for models with a degenerate ground state space. Hence, 
at this point, we can certify the overlap between the true 
state and a small subspace of the exponentially large am- 
bient space. 

It remains to be shown that even this comparatively 
small ambiguity can be efficiently resolved. Treating all 
special cases which may appear in the most general sit- 
uation is somewhat cumbersome and will be deferred to 
a future publication. Here, we restrict attention to the 
most relevant case: a two-fold degeneracy (this covers 
the GHZ and the W state). Formally, we are facing the 
task of performing tomography in a copy of C 2 , which 
has been embedded into C d . A classic result states that 
any two orthogonal pure states of a multi-partite sys- 
tem can be reliably distinguished using local operations 
alone [17]. This holds true in particular for the embedded 
versions of the eigenvectors of the three Pauli matrices 
acting on C 2 . Therefore, the overlaps between a state in 
the ground-space and this informationally complete set 
of vectors may be found by local operations. We end the 
discussion by noting that the procedure equally applies 
to the case where p is far from a pure state, as long as 
its support is contained in an MPS space. 

Generalizations and Outlook. The theory of MPS 
can be generalized to higher-dimensional configurations, 
where the resulting states are sometimes referred to as 
PEPS. We note that our gap estimate does not cru- 
cially rely on the geometry of the neighborhood relations. 
Therefore, conceptually, it should be relatively straight- 
forward to generalize the basic techniques to these situ- 
ations. However, efficiently finding good PEPS approxi- 
mations to the data is likely to be extremely challenging. 

While this Letter emphasized (approximately) pure 
states, tomography of highly mixed states may actually 
be more pertinent in realistic settings. If we drop the de- 
sire to certify the result, this task proves relatively simple. 
Recall that a mixed state p is called finitely correlated [2] 
if the linear space of operators 



{tr R { P A R ]\A R eB(R)} 



(7) 



has dimension smaller than some fixed constant D 2 , 
where the indices L, R refer to a bi-section of the chain 
into a left and a right part. The symbol B(R) denotes 
the set of all operators acting on the right hand side 
(r.h.s). Physically, the definition says that there are only 
D 2 ways of modifying the state on the l.h.s. by condition- 
ing on an event taking place on the r.h.s. of the chain. 
It has been shown in [2] that any state with that prop- 
erty has an efficient representation, and once again it is 
natural to ask whether this small set of parameters can 
be experimentally obtained, and whether success can be 
certified without technical assumptions. Under a mild in- 
vertibility assumption akin to (4), the answer to the first 
question turns out to be affirmative. While details will 
be provided elsewhere, we briefly list the main features of 
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the algorithm: It is guaranteed to obtain a faithful repre- 
sentation, without the need for any optimizations or vari- 
ational calculations; it yields a representation which may 
be used to predict the expectation value of any (of the 
exponentially many) factorizing observables, even though 
only (linearly many) local observables were measured to 
find the estimate. The MPS pure states discussed above 
are included as a special case. On the downside, the 
resulting "matrix product operators" are much harder 
to interpret; there is no way of heralding errors if the 
assumptions are not met and the resulting state will be 
non-negative only within the margins of experimental ac- 
curacy. 

We remark that our the results on mixed state tomog- 
raphy immediately generalize to quantum process tomog- 
raphy. The analogue of a non-factoring quantum state 
is here a channel with correlations over many uses, while 
the analogue of the FCS representation are channels with 
memory [18]. Obtaining the parameters of these objects 
efficiently is, courtesy of the Choi-Jamiolkowski isomor- 
phism, just as straight-forward as learning mixed FCS. 
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Note added. After the research for this work had been 
completed [19], the authors learned that two independent 
groups had worked along similar lines. 

Let us compare the present work with the recent work 
in Ref. [20]. Conceptually, the two contributions seem 
very close. Both make the point that tomography of 
many-body quantum systems may not be as impossible 
as it superficially seems. But the respective technical fo- 
cus is quite different. While Ref. [20] puts forward a new 
algorithm for finding MPS approximations given mea- 
sured data, we emphasize methods for certifying that the 
estimate is actually correct. Employing the MPS-SVT 
method of Ref. [20] in step one of our H-MPS algorithm, 
it should be possible to seamlessly combine the two re- 
sults. 

Unrelatedly, the authors were recently made aware of 
the fact that a further group was preparing results [21] 
on the problem presented here. 



